Mixing-induced global modes in open active flow 
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We describe how local mixing transforms a convectively unstable active field in an open flow into 
absolutely unstable. Presenting the mixing region as one with a locally enhanced effective diffusion 
allows us to find the linear transition point to an unstable global mode analytically. We derive the 
critical exponent that characterizes weakly nonlinear regimes beyond the instability threshold and 
compare it with numerical simulations of a full two-dimensional fiow problem. 
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In many natural and laboratory flows active chemi- 
cal and biological processes may occur. Examples in- 
clude chemical reactions in minimixers, plankton growth 
in ocean, see e.g. [ij and references therein. Quite often 
this activity occurs in an open rather than in a closed 
geometry. Here the main issue is whether the through- 
flow is stronger or weaker than the activity. One has to 
compare the velocity of the throughflow with the velocity 
of the activity spreading due to diffusion. If the through- 
flow is stronger, the activity is "blown away" like a candle 
flame in a strong wind, in the opposite case a sustained 
activity can be observed [3| ■ This simple picture is valid, 
however, only for homogeneous media. Often additional 
vortexes are superimposed on a constant throughflow, 
due to, e.g., wakes behind islands in ocean currents or 
mixing enforced by revolving fan blades in laboratory 
experiments. We want to study under which conditions 
such an additional kinematic mixing in a strong open 
flow can lead to a transition to a sustained activity, and 
to characterize this transition quantitatively. 

Our main model is a reaction-advection-diffusion equa- 
tion for the dimensionless concentration of an active 
scalar fleld u{r, t) 

Mf + (V + W(r, t)) ■ Vu = DoV^u + au{l - u^) . (1) 

Here V = (V, 0, 0) is a constant throughflow in x- 
direction, Dq is a molecular diffusion of the scalar field. 
Activity is assumed to be of the simplest form: a linear 
growth with rate a with a saturation at u = 1. The non- 
linearity index p is typically integer (1 or 2) for chemical 
reactions, while for biological populations a wide range 
of values of p has been recently reported Q . Mixing is 
described by a spatially localized incompressible veloc- 
ity field W(r, i), its intensity is denoted as W. Note 
that in the absence of fluid flow Eq. ^ is reduced to the 
famous Kolmogorov-Petrovsky-Piskunov-Fisher (KPPF) 
model of an active medium with diffusion (see, e.g., 4] for 
original references, analysis, and applications of KPPF), 



while for a = Eq. ([T]) describes a linear evolution of a 
passive scalar in a flow. Model ([T]) can be used for the 
description of biological activity, where u is, e.g., concen- 
tration of a growing plankton, advected by oceanic cur- 
rents; for a possible laboratory realization see recent ex- 
periments with an autocatalytic reaction in a Hele-Shaw 
cell with a throughflow p| . 

In the absence of flow, the diffusion causes the active 
state to spread forming eventually a front with velocity 
Vf = 2\fal\ Thus, for vanishing mixing = 0, 
the activity is blown away provided V > Vf. For this 
parameter range the instability in Eq. ([1]) is convective 
and in the absence of external sources, no activity is ob- 
served in the medium. A nontrivial state is, however, 
possible if there is a localized source of the field u: then 
a growing tail stretches from this source in the down- 
stream direction, where it eventually saturates at u = 1. 
The linearized problem with a point ((^-function) source 
of intensity s can be readily solved, yielding 

e{V)-^ exp[xV/2Do] exp{~\x\V /2Dq) 

uir) = ei2TrDo)~' eMxV/2Do]Koi\r\V /2Do) 



* E-mail: arthur.straube@gmail.com 
Paper published in Phys. Rev. Lett 99, 184503 (2007) 



(2) 
(3) 

in one- and two-dimensional setups (where Kq is the 
modified Bessel function, V ^ (V'^ - Vf)^^^). Note that 
in both solutions u ^ exp(fiz^x) a,s x ±oo, where 
M± = i2Do)-HV±V). 

In this letter we demonstrate, that beyond some crit- 
ical intensity Wcr, a localized mixing W(r,i) turns the 
convective instability locally into the absolute one, which 
results in a stationary (in statistical sense) proflle of u 
(see Fig. [1] for a sketch of the proflle and Fig. |4] below 
for a numerical example). Beyond criticality W > Wcr, 
the mixing region acts as an effective source of the field, 
in Fig. [1] this region is denoted as a "source." A "tail" 
where the field grows exponentially as in ([2]), ([3]) extends 
downstream of the source. Further downstream stretches 
the "plateau" domain where u = 1. Our main quantita- 
tive result, obtained by matching solutions in these three 
domains, is the critical exponent /? relating the intensity 
of the effective source £e// with the mixing intensity: 

Seff^iW-Wcr^. 

To develop the theory we use the concept of global 
modes |7i]. In this concept a self-sustained non-advected 
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(a) "source" "tail" "plateau" 



situation, described by the linearized at w = Eq. ([4]): 




FIG. 1: (a): Quasi-one-dimensional open flowwith a localized 
mixing zone. Panels (b) and (c) illustrate the construction of 
the nonlinear global mode, they show qualitative profiles u{x) 
and z{x) = ^ in the linear approximation at the criti- 
cality (dashed line) and in nonlinear regime slightly beyond 
criticality (full line). In the latter case the profile is nearly 
linear for x < Xr but deviates due to nonlinear terms for 
X > Xr, see discussion of Eq. (I14|l . Regions "source," "tail," 
and "plateau" are explained in the text. 



pattern arises due to inhomogeneities of the system. Typ- 
ically, one considers inhomogeneities of the growth rate 
a: if a = a(r) has a hump where locally the front ve- 
locity is large Vf loc > V, then a global mode appears, 
located at this hump and downstream. In this paper 
we are interested in another, mixing-based mechanism 
of a global mode appearance. It can be easily under- 
stood if the concept of effective diffusion (see, e.g., Q) 
is used to describe the mixing term in ([l]). In this ap- 
proach, we phenomcnologically introduce effective diffu- 
sivity D{r) — Dq + D„iix(r) that accounts for the coarse 
grained mixing dynamics, and write instead of ([T]) an 
equation with a non-homogeneous diffusion 



ut + V ■Wu = V[L»(r)Vu] + au{l - u^) 



(4) 



A hump of diffusivity D{r) leads to an increase of the lo- 
cal front velocity Vf , and one expects that when the front 
propagation prevails over the throughflow, a stationary 
global mode can appear, producing a mixing-induced 
sustained structure. We focus on a geometry shown in 
Fig. [IJa): in a constant open flow there is a localized 
region of strong mixing, which, as we will see below, is 
not necessarily chaotic or turbulent. The theory below 
will be developed for a one-dimensional case, which is 
relevant, e.g., for flows in a micropipe; the results will be 
supported by numerical simulations of two-dimensional 
flows. We restrict ourselves to this case because of com- 
putational simplicity, and also because two-dimensional 
flows are relevant for many geophysical and laboratory 
(especially in microfluidics) experimental situations. 
We start with a linear analysis of a one-dimensional 



dt dx dx 



(5) 



We look for an exponentially growing in time solution 
and with an ansatz u{x,t) ^ exp[At -I- / z{S,)d^] obtain 



dz 

dx 



V - 



dPjx) 
dx 



D{x) 



D{x) 



(6) 



As 
D 



x\ oo we have a homogeneous medium with 
Dq, here the solution should tend to values = 

at which the r.h.s 



of 





{2Do)-HV ± y/V^ - 4(a - X)Do 

([6|) vanishes. More precisely, as a; — > —oo we have z z 
and as a; — > oo we have z — s- With these two bound- 
ary conditions one easily finds the solution of ([6]) numer- 
ically, matching at z = integrations starting at large 
\x\ from the values z^.. In a particular analytically solv- 
able case of a piecewise-constant diffusivity: D = Dq for 
|a:| > I and D = Di > Dq for \x\ < I, one can perform 
the integration analytically and obtain the equation for 
the growth rate A: 



I = 



2Di 



arctan / 



I V^ - Aja- X)Do 
Aia~X)Di - 



(7) 

The value A = corresponds to the onset of global in- 
stability, in this case ([7]) gives the relation between the 
critical values Icr and Dicr- 



2D^c 



: arctan - 



- Vf 



1/2 



(8) 



From ([5]) it follows that Dicr — > OO as Icr Imin = 

V/{2a). In other words, there exists a minimal size of 
the mixing region, so that for I < Imin even a very strong 
mixing, with a very large effective diffusion, cannot create 
a global mode (the same is true for a smooth profile of 
D{x); note also that the size of the mixing region is not 
limited from above). This is in contrast to the situation 
when the global mode is induced by a local hump of the 
growth rate a (cf. Q): here one can obtain instability 
even when a{x) is highly localized (a delta- function) , a 
global mode then looks as in ((21), A similar analysis 
performed for a two-dimensional inhomogeneous version 
of Eq. ([5]) also yields a minimal radius of a diffusive spot 
that can lead to instability. 

Next we discuss the linear stability not in the frame- 
work of the effective diffusion model but in the full 
reaction-advection-diffusion problem as in Eq. ([T]) . After 
the linearization we arrive at a linear stability problem 



Ut + (V + W(r, t)) ■ Vu ^ DoV^u - 



au 



(9) 



which is nonstationary if the velocity field W is time 
dependent. Then the proper way to determine the sta- 
bility is to calculate the largest Lyapunov exponent (LE) 
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A = In I |u| 1). This can be done numerically, as de- 
scribed in ref. [10| . Noteworthy, in this consideration we 
are not restricted to a deterministic flow, as the LE can 
be calculated also for a randomly or chaotically time- 
dependent field W(r, i). 

We have calculated the LE for a linearized two- 
dimensional reaction-advection-diffusion equation ([9]) 
subject to a constant open flow with a superimposed os- 
cillating vortex, described by the stream function 

'ifix,y) ^ Vy + W exp[~{x^ + y^)R-^]cos{Ljt) . (10) 

We fixed V = 1, R = 1 and cu — 2, and calculated the 
LE A for different molecular diffusion constants Dq and 
vortex intensities W (Fig. [J)- Note that the parameter a 
simply shifts the LE, therefore we plot A— Aq, where Aq = 
a — V'^/{4Dq) is the LE for a non-mixed flow. One can 
see that the mixing-induced enhancement of field growth 
is mostly pronounced for small diffusion and is maximal 
at « 3. This is the mixing strength at which a chaotic 
saddle [l[ in the Lagrangean particle trajectories appears. 
A further increase of the vortex intensity does not lead, 
however, to a significant growth of the LE. 

Now we develop a nonlinear theory of the global mode. 
It is clear that the nonlinear saturation stops the expo- 
nential growth of a slightly supercritical linear mode and 
leads to a nonlinear solution with a finite amplitude. Our 
aim is to describe the dependence of this amplitude on 
the deviation from criticality. First we notice, that the 
very notion of the amplitude is here nontrivial. Indeed, 
the nonlinear solution looks as in Fig. [TJj (cf. Fig. [3] be- 
low); it saturates to u = 1 in the downstream direction. 
However, outside the mixing region the field looks like a 
solution caused by a localized field source. Thus, we can 
take the effective intensity of this source Seff, which is 
proportional to the characteristic field amplitude in the 
mixing region u(0) [see relations ©, Q], as the order 
parameter of the transition. The deviation from the crit- 
icality we will measure with the growth rate A, for which 
holds X (xW — Wcr in the full model ^ or X oc D — Dcr 
in model 

We will consider the simplest possible setup, namely 
the nonlinear modification of one-dimensional Eq. ([5]): 



du 
'dt 



V 



du 
dx 



d_ 

dx 



D{x) 



du 
dx 



au{l 



(11) 



We look for a stationary global mode m(x), and rewrite 
this equation as the system 



dz 

dx 
du 

dx 



V 



dPjx) 
dx 



D 



a 
D 



D 



(12) 
(13) 



We consider this system separately in two spatial do- 
mains. The first, linear region, includes the inflow and 
the mixing domains ("source" in Fig.[T]): —oo<x<Xr, 




FIG. 2: (Color online) Lyapunov exponent characterizing sta- 
bility of the global mode mixed by periodically blinking vortex 
(jlOp . A similar picture holds for a stationary vortex. 

where the field u{x) remains small. In the second, out- 
flow region Xr < x < oo, the field u further grows ("tail") 
and nonlinear ly saturates ("plateau"). In the linear re- 
gion, because of smallness of the field, we can neglect 
in (fT2|) . thus we obtain an equation similar to The 
only difference is that because we look for a stationary 
solution, in IT^ the term ~ A is absent. Near the crit- 
icality, where A is small, we can consider this term as a 
perturbation, therefore the solution of (|12[) in the linear 
region is close to the solution of Eq. it has the asymp- 
totic z as x — > — CX3. Due to the perturbation term 
(X A, at the right border of the linear region z deviates 
from the deviation /i_ — z{xr) is proportional to A, 
and, thus, to I? — Dcr- At Xr the field u is small and 

u{Xr) fx u{0). 

Next we consider full equations p^ . in the non- 
linear region x > x,.. Here the solution should tend at 
X — > oo to the saddle fixed point u — I, z — 0. Thus, 
starting integration from large values of x in the nega- 
tive direction, we have to follow the stable manifold of 
this saddle and match this solution dA, x — Xr with the 
obtained above. Because the value to be matched z(xr) 
is very close in the region where the solution (z,u) 



approaches {z{xr), u(xr)) we can write 
to obtain 



dx 



Az = (/i_|_ — /i_)Az - 



D 



uP{xr)e- 



(14) 



Here, since Az = /i_ — z(x) is small, we have approxi- 
mated the solution of as u ~ u{xr)e'^-^^~'^^\ Be- 
cause linear inhomogeneous Eq. (jl4p is solved in the neg- 
ative in X direction, the solution follows the slowest expo- 
nent: Az oc ex.p[j{x—Xr)], where 7 = min(/i+— /i_,p/i_). 

At the criticality, the region of validity of the ex- 
ponential solution Az cx exp[7(a; — Xr)] becomes very 
large. Thus it is dominant for small deviations from crit- 
icality A, therefore we can estimate the coordinate Xg 
at which the field u saturates (i.e., we reach the state 
u 1 and z « 0) from the relations above: from 
— ^_ w {z{xr) — ^-)e'^^^'~^''^ it follows {xs — Xr) ~ 
— 7^^1n(/i_ — z{xr))- Substituting this in the expres- 



sion for u{x), we obtain u{xr) oc (/i_ 
we take into account that /i_ — z{xr 



-z{Xr)r-^^. 
OC D ~ Drr- 



Now 
and, 
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FIG. 3: (Color online) The critical exponent /? calculated for 
the model ^ (symbols), compared with theoretical predic- 
tion (|16|) (lines). One can clearly see the crossover between 
two regimes of the field saturation in dependence on param- 
eter a, the latter is related to v in (|16p via v ~ a~^^^ . 

because the evolution of u in the interval Q < x < Xr only 
weakly depends on the criticality, Seff = w(0) cx u{xr)- 
The final expression for the scaling law of the amplitude 
of the global mode thus reads 

Eeff OC A'', /3 = — = max( — — ,-]. (15) 

The critical index [3 depends only on the nonlinearity 
index p and on the dimcnsionlcss velocity v = V/Vf: 

I iif^S if 1< w < . ^ ^ 

This main result of our letter can be physically inter- 
preted as follows. The exponent /3 is determined solely 
by the nonlinearity index p if the throughflow velocity 
is much larger than the front velocity {v large). Here 
the field in the plateau domain (see Fig. [T]) is effectively 
uncoupled from the source, and the saturation of the in- 
stability is due to the local nonlinearity at the source. For 
a small throughflow velocity {v close to one) the plateau 
state interacts with the source via the tail. Due to this 
"remote control," the field at the source is saturated more 
efficiently than due to nonlinearity, here the exponent f3 is 
determined solely by the form of the tail, which depends 
on the velocities ratio v. 

Below we check formula (I16p with direct numerical 
simulations of model ([T]). A stationary vortex (jTU]) with 
Lu = and R = 1 was imposed on a constant flow with 
V ~ 1. Keeping the diffusion constant flxed Dq — 0.3, 
for different field growth rates a we have found, from the 
linearized equations, the critical vortex intensities Wcr 
at which the global mode becomes first unstable. Then 
we solved full nonlinear equations close to criticality and 
found the exponent /? according to p3|) . The stationary 
problem was solved with a finite difference method in a 
domain 0<a;<60,0<y<40 with periodic boundary 
conditions in y and conditions u(0) = 0, §f (60) = 0. The 
results are presented in Fig. [3l they are in a very good 
agreement with the theoretical prediction (fTS]) , (fTB]) . Fig- 
ure m shows the example of the stationary mode appear- 
ing beyond the instability threshold. 




FIG. 4: (Color online) The active field behind a vortex with 
W ^ 4 placed at a; = y = 0, for 1/ = 1, a = 0.5, Do = 0.3, 
p = 2. 



In summary, we have described the mixing-induced 
transition from a convectively unstable active field in 
an open flow to a persistent global mode. Our theo- 
retical approach bases on the representation of the mix- 
ing region as a domain with locally enhanced effective 
diffusion. For a nonlinear regime close to criticality we 
have derived the critical exponent /? that depends 
only on two parameters of the system: the dimension- 
less flow velocity v normalized by that of the front, and 
the nonlinearity index p. For large velocities the criti- 
cal exponent depends only on the system's nonlinearity, 
which means a local in space saturation of the instabil- 
ity. For small velocities the exponent is a function of 
velocity, here the growing downstream tail of the active 
field imposes the saturation. Notably, this prediction of 
the one-dimensional theory is in a good accordance with 
two-dimensional calculations. 
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